Three-particle integrals with the Bessel functions. 

Alexei M. FroloJ] 

Department of Chemistry 
University of Western Ontario, 
London, Ontario N6H 5B7, Canada 
(Dated: January 15, 2013) 

Abstract 

Analytical formulas for some useful three-particles integrals are derived. Many of these integrals 
include Bessel and/or trigonometric functions of one and two interparticle (relative) coordinate 
r 32) r 3i and ?*2i. The formulas obtained in such an analysis allow us to consider three-particle 
integrals of more complicated functions of relative/perimetric coordinates. In many actual prob- 
lems such three-particle integrals can be found in matrix elements of the Hamiltonian and other 
operators. 
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I. INTRODUCTION 



The main goal of this study is to consider some special integrals which are closely related 
with various fundamental three- and few-body problems in physics. The are of possible 
applications of such integrals include atomic, molecular, nuclear and other areas of physics. 
In many of the problems from these areas one finds similar integrals which must be taken 
over three scalar distances and these three distances correspond to the sides (or ribs) of the 
triangle formed by three 'particles'. On the other hand, it is clear now that the general 
theory of three-particle integrals is a rapidly growing area of mathematical physics. The 
methods developed for numerical evaluations of different three-particles integrals can be 
used in various mathematical problems. 

In general, the three-particle (or three-body) integral designates the following integral 

I(a, (3, 7; F) = J J J F(r 32l r 31 ,r 2 i) exp(-ar 32 - (3r 31 - ir 21 )r 32 r 31 r 21 dr 32 dr 3l dr 2l (1) 

where a, f3 and 7 are the three real values which are usually called and considered as the 
non-linear parameters. The function F(x,y,z) in Eq.flT]) is the analytical function of each 
of the three real variables x,y,z. The generalization of Eq.flT]) to the case of the complex 
variables is possible, but in this study we do not want to discuss it. Moreover, to simplify 
the problem below we shall assume that the function F(r 32 , r 3 %, r 2 \) does not grow 'very 
rapidly' when its arguments increases to the infinity. Such an assumption is needed to 
guarantee the convergence of all arising integrals (see below). The variables in Eq.flT]) 
are the three relative coordinates = | r j — Vj | . These coordinates are the scalars and they 
are symmetric upon permutation of their indexes r\y = r^. The relative coordinates are 
not truly independent of each other, since for these coordinates the six following conditions 
must always be obeyed: >| — | and < + rjk, where k) = (1,2,3). 

The integrals, Eq.flT]), can be found in various few-body problems, but in the Quantum 
three-body problems they play a central role, since the expressions for all matrix elements 
of the Hamiltonian and overlap matrices are reduced to Eq.flT]). Furthermore, the explicit 
expressions for all expectation values (regular and singular) also reduced to the computation 
of formulas each of which coincides with Eq.flT]). Different approaches to analytical and 
numerical computation of the three-body integrals Eq.flT]) were developed in the middle of 
1980's [lj], Q]. In particular, in those works it was shown that the most convenient and 
simple way to compute such integrals is based on the use of the three perimetric coordinates 
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Ui,u 2 ,u 3 |3J, |4| instead of relative coordinates r 3 2,r 3 i,r 2 i mentioned above. The relation 
between these two set of coordinates is linear: 

ui = -(r- 3 i +r 2 i -r 32 ) , w 2 = -(r 21 + r 32 - r 31 ) , it 3 = ^(r 32 + r 3i - r 21 ) (2) 

The inverse relation takes the form rij = Ui + iij. The use of the three perimetric coordi- 
nates have a number of advantages in calculations of various three-particle integrals. For 
instance, the three perimetric coordinates Ui,u 2 ,u 3 are independent of each other and each 
of them varies between and +00. The substitution (r 32 , r 3 i, r 2i ) — )■ (u\,u 2 ,u 3 ) drastically 
simplifies analytical and numerical computations of all three particle integrals. In perimetric 
coordinates ui,u 2 ,u 3 the basic integral, Eq.([T]), is written in the form 



exp[-(a + (3)u 3 - (a + ^)u 2 — (/3 + 7)^1] (wi + u 2 )(ui + u 3 )(u 2 + u 3 )du 1 du 2 du 3 (3) 

where the factor 2 in the front of the integral is the Jacobian of the (r 32 , r 3 i, r 2 i) — > (u\,u 2 , u 3 ) 
transformation. 

The integral, Eq.(j3]), takes the form of the Laplace transformation. It is also clear that 
the image function I will be written in the form I(a + j3, a + 7, (3 + 7;^). In general, 
the calculation of the three-body integrals, Eq.Q, appears to be closely connected with 
the Laplace transforms. Indeed, the tables of Laplace transormations are of great help 
in analytical and numerical computations of various three-body integrals. In our earlier 
studies we have derived the explicit formulas for a large number of three-body integrals. 
Such formulas include different regular and singular integrals, integrals with logarithmic 
terms, etc. However, the formulas for some important three particle integrals have never 
been derived in earlier studies. For instance, the integrals which contain one or two Bessel 
function, B and inte^ ta which one function of Native coc,d,„ at e i 8 .presented in SOm e 
'difficult' forms, e.g., as the infinite power series, or as approximate expansions written in 
terms of other functions. In this work we consider some of such integrals. 

II. EXPONENTIAL VARIATIONAL EXPANSION 

As mentioned above the three-particle integrals, Eq.(TjQ), arise in various three-body prob- 
lems. In general, all problems related with the construction of highly accurate approxi- 
mations of the actual wave functions for bound states lead to such integrals. This means 



I(a, /3, 7; F) = 2 / / / F{u 2 + u 3 ,ui + u 3 ,ui + u 2 ) 




x 
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that the three-particle integrals, Eq.flT]), can be found in all bound state problems, including 
transitions between bound states, time-evolution of bound states, etc. On the other hand, 
analytical/numerical computation of the integrals, Eq.([T]), is the central part of problems 
related with the photodetachment and decays of bound states in three-body systems. Fortu- 
nately, the three-particle integrals, Eq.([T]), can be considered as the matrix elements between 
two exponential basis functions written in the relative coordinates. In general, the exponen- 
tial variational expansion of the three-body wave function for the bound S(L = 0)— states 
is written in the form 

1 N 

V = -(1 + kP 2 i) °i exp(-air 32 - &r 31 - 7^21) (4) 
1 i=i 

where k = ±1 for symmetric systems and k = otherwise (see below). The generalization 
of Eq. (HJ) to the case of bound states with arbitrary angular momentum L is written in the 
form [6] 

1 N 

* = r(l + kP 21 ) £ £ Cy^(r 31 , r 32 ) exp(-«,r 32 - &r 31 - 7i r 21 ) (5) 

1 1=1 lx 

jut much more flexible!) generalization of Eq.(j4]) to the case 



A slightly more complicated 
of arbitrary L takes the form 

1 N 

^lm = -(1 + «Ai) Y.Y. C iyLM 2 ( r 3i,r32)M r 32,r 31 ,r 2 i)exp(-a i u 1 - (5iU 2 - 7^3) x (6) 

1 i=l lx 

exp(z5jMi + ie { u 2 + ifiU 3 ) , 

where Cj are the linear (or variational) parameters, aij, 7$, Si, ej and /j are the real non- 
linear parameters and z is the imaginary unit. In the last equations all exponents contain 
the three perimetric coordinates ui,u 2 , u 3 instead of relative coordinates r 32 , r 31 , r 2 i used in 
Eq.®. The functions 3^ 2 (r 3 i,r 32 ) in Eqs.flH]) and © are the bipolar harmonics |8|] of the 
wo vectors r 3 i = r 3 i • n 3 i and r 32 = r 32 ■ n 32 . The bipolar harmonics are defined as follows 



3^ 2 (x, y ) = xV 2 E ^Zx-,hrn 2 Y hmi (n x )Y £2m2 (n y ) (7) 

where C^[^ <f m are the Clebsch-Gordan coefficients (see, e.g., [8|) and the vectors n x . = - 
and n y = ^ are the corresponding unit vectors constructed for arbitrary non-zero vectors x 
and y. Also, in this equation L is the total angular momentum of the three-body system, i.e. 
L 2 ^ lm = L(L + l)ty lm, while M is the eigenvalue of the L z operator, i.e. L z ^ L m = M^lm- 



In actual calculations it is possible to use only those bipolar harmonics for which l\ + l 2 = 
L + e, where e = or 1. The first choice of e (i.e. e = 0) corresponds to the natural spatial 
parity xp = (~ 1) L °f the wave functions jj. The second choice (i.e. e = 1) represents states 
with the unnatural spatial parity xp = ( — 1) L+1 - In almost all works on highly accurate 
bound state computations only the bound states of natural parity are considered. In real 
physical systems only such states exist and stable. 

The polynomial-type functions </>i(r 32 , r 3 i, r 2 i) are used in Eq.flT]) to represent the inter- 
particle correlations at short distances. In general, such simple polynomial functions al- 
low one to increase the overall flexibility of the variational expansion Eq.flT]). In many 
studies, however, these additional functions are chosen in the form 0i(r 32 , r 31 , r 21 ) = 1 for 
i = 1, . . . , N, since the overall convergence rate of the variational expansion, Eq.flH]), is al- 
ready very high. The operator P 2 i i n Eq.flT]) is the permutation of the identical particles in 
symmetric three-body systems, where k = ±1, otherwise k = 0. 

Note that the basis functions, Eq.fl!]), are the (2L+l)-dimensional vectors, while all matrix 
elements of the Coulomb three-body Hamiltonian matrix and overlap matrix are scalars. In 
reality, one finds no contradiction here, since the angular integral of the products of the two 
and three bipolar harmonics 3^i,M 2 ( r 3i> r 32) equals to the products of some linear functions of 



the r| 1; r| 2 , r 3i • r 32 variables and 6j— and 9 j— symbols, respectively (see, e.g., [lOj and [11]). 
The scalar variables r| 1 (= u 3 + m^ + 2m 1 m 3 ), r\ 2 {= 143 + ^2 + 2^%), r 31 ■ r 32 = § (^"32 + r li — r ii ) 
are easily expressed as the quadratic functions of the perimetric coordinates. Indeed, by 
using the relations between the relative and perimetric coordinates one finds 

r 3i = u l + u i + 2wi^ 3 , r\ 2 = u\ + u\ + 2u 2 us 
r 3i ■ r 32 = i(rj 2 + r\ x - r^) = u\ + u x u z + u 2 u 3 - u x u 2 (8) 

Numerical computations of the corresponding radial integrals are slightly more complicated 
than in the case of L = 0. However, all integrals needed in actual bound state computations 
with the use of the exponential variational expansion, Eq.flE]), are written in the form of 
Eq.flTJ, or can be reduced to such a from. This fact explains our permanent interest in 
developing of the new analytical/numerical approaches to calculations of the three-body 
integrals, Eq.flTJ). 
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III. INTEGRALS OF POLYNOMIAL FUNCTIONS 



First, let us present here our formula for the three-body integral, Eq.flTJ), which includes 
the polynomial function F(r 32 , r 3 i, ^21) — r 32 r 3i r 2i- I n this case the arising integral is desig- 
nated as r fc; ; ;n (a, j3, 7) and it is written in the form 



I\z ; „(a, ft, i) = J J J r 32 r 31 r r ^ exp(-ar 32 - 0r 31 - ^r 21 )dr 32 dr 3l dr 2l (9) 

where all indexes k, I, n are assumed to be non-negative. In perimetric coordinates the 
integral, Eq. ([§]), is written in the form 

r+oo r+oo r+oo 

r fc;Z;n (a,/3, 7 ) = 2 / / / (u 2 + u 3 ) k ( Ul +u 3 ) l ( Ul +u 2 ) n x 
Jo Jo Jo 

exp[— (a + /3)u 3 — (a + ^y)u 2 — (f3 + j)ui]duidu 2 du 3 (10) 



Analytical calculations of the integral, Eq. fllQI) . is simple. Finally, one finds the following 
formula 

r f,f),i = 9Vvf r'r'r (l-h + hy. (k-h + n^l (n-m + ZQ! 

k;lA ' P,1) k % to nh 1 h 1 (« + P) l - h+k ^ +1 (a + 7 )*-*i+»i+i (0 + 7 )«-»r«i+i 

fc l n s~<k\ (~<ri\ 

= 2 A-l /' r;' ^n-m+k^k-lg+h^l-h+m /-.-.x 

' ' '^ ^on^o(« + ^)'~' 1+/ci+1 (a + 7) /c " fcl+ril+1 (^ + 7)^" 1+ ' 1+1 

where C^} is the number of combinations from M by m (m and M are the integer non- 
negative numbers). The formula, Eq. ffTT]) . can also be written in a few other equivalent 
forms. The function x ^ +1 in Eq. ffTT]) is the A n (X) function introduced by Larson [12|. For 
the first time, I have produced the formula, Eq. ffTT]) . in the middle of 1980's (see [Jl^] and 
references therein). 

Note that in some our earlier works the following integral in perimetric coordinates was 
considered as the basic three-body integral: 

. . r +co f+°° f +0 ° u p l 1 u 2 2 u p ' i exp(-a-ui - bu 2 - cu 3 ) . 

B(a,b,c;p 1 ,p 2 ,p 3 ;q ,q 1 ,q 2 ,q 3 ;s) = / / ■ ■ ■ du ± x 

Jo Jo Jo (g + QiUi + q 2 u 2 + q 3 u 3 ) s 

r(pi + l)r(p2 + l)T(p 3 + 1) r°° exp(-q x)x s - 1 dx 

du 2 dU 3 = —7- -r / 7 r — r — r — ,12 

r(s) Jo (a + q 1 x)Pi+ 1 (b + q 2 x)Pi +1 (c + q 3 x)Pt +1 V ; 



where T(x) is the Euler's gamma-function (see, e.g., |l3j, 14j). This integral depends upon 
eleven parameters a, b, c,p%,p 2 ,p 3 ] q , qi, q 2 , q 3 and s. Formally, all these parameters must be 
real and positive (or non-negative). In particular, for s = 1, go — 1 and gi = q 2 = q 3 = one 
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finds from the last formula 

B(a, b, c; Pl ,p 2 ,p 3 ; 1, 0, 0, 0; 1) = aPl+1&p2+ltf , 3+1 (13) 

This formula leads to the generalization of Eq. ljTTj) to the case of non-integer values of k, I, n. 
It is often used to operate with the modified basis sets, e.g., with the basis set which includes 
semi-integer powers of perimetric coordinates. 

In some three-body problems related, e.g., with the scattering, one finds a number of 
advantages in using of some non-exponential basis sets, e.g., power- type wave functions of 
the relative and/or perimetric coordinates. In such cases we need to determine different 
basic integrals. In this study we do not want to discuss this interesting problem. Note only 
the following formula which arises in the case of of power-type basis functions 

,+oo r+oo r+oo ul l u p 2 2 u Pi 

(g + qiu\ + q 2 u 2 + q 3 u 3 )' 



,+oo ,+oo ,+oo uVvZ'vr* , , , 

G{Pi,P2,P3; qo,qi,q2, s) = / / / - — ■ ■ ■ -du 1 du 2 du 3 

Jo Jo Jo (g + qiui + qiu 2 + g3^3j 

r(pi + i)i> 2 + i)r(p 3 + 1) r(a - Px - P2 - p 3 - 3) 



(14) 



F('c.^ Pl + 1 /YP 2 + 1 /rP 3 + 1 S-P1-P2-P3~3 

1 l s Jyi q 2 q?, yo 
where it is assumed that s > p\ + p 2 + p 3 + 3 and all values pi, qi (i = (1,2,3)) and go must 
be positive. 



IV. DERIVATION OF THE RELATED INTEGRALS 

By using the expression for the I\z ;n integral, Eq.(|9]), we can obtain analytical formulas for 
various three-particle integrals. First, consider the matrix elements of the real (analytical) 
functions f{r 2 \) and F(r 2 \) which are represented by the following series: 

f(r 32 ) = J2 A n r 32 and F(r 32 ) = J2 A nr% 2 exp(-B n r 32 ) (15) 

n n 

where the number of terms can be finite, or infinite. The computation of the matrix elements 
of these functions is reduced to the calculation of the two following sums 

M / = ]TA n r n+1;1;1 (a,/3,7) (16) 

n 

and 

M F = J2 A n V n+1;1 . A (a + B n , f3, 7 ) (17) 
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where the integals T n -k-i are defined in Eq.flT]). If the coefficients A n in these expansions 
rapidly decrease with n, then one needs to compute only the finite number of terms in 
such sums. For instance, if the coefficients in Eq. (ll5j) are the power-type functions of some 
small parameter, then the series in Eq. ffTET) and Eq. ( I17j) converge rapidly. This is the case in 
various atomic problems related with QED applications, when A n ~ a n , where a = j- ~ ^ 
is the dimensionless fine-structure constant. 

This approach can also be used to produce the analytical formulas for more complicated 
integrals, e.g., the general three-body integrals with the Bessel function j^Vr^)- First, let 
us obtain the computational formula for the following integral 

5 S;n(«> P> T, V) = J J J r 32 r 3i r 2iio(^ 32 ) exp(-ar 32 - (3r 31 - ~fr 21 )dr 32 dr 31 dr 21 (18) 
By using the formula jo(x) = and Eq. ffTo 7 ]) one finds 

oo (_l\qy2q q max ( _l)qy2q 

*2«(«, 0, 7; V) = £ \ 2q+1 y T W^ a > P> T) ~ E \ 2q + 1} , r Wn(«, 0, 7) (19) 

The integral B^J. n (a, (3, 7; V) in the last equation converges for all V, but for V < 1 it 
converges very rapidly. In reality, the maximal value of the index q (or q max ) in Eq. ffl9|) is 
finite. Numerical investigations indicate that to stabilize 15 decimal digits for V < 1 one 
needs to use in Eq.f JT9|) q max = 20 — 40. For V > 2 the value of q max rapidly increases up to 
50 - 70 and even 100. The same conclusion is true about the convergence of the three-body 
integrals with the Bessel function ji(x) = ~ This integral takes the form 

Bll n (a, (3, T , V) = f ^^|j^ r fc+2g+1;;;rt («, {}, 7 ) (20) 

The three-body integrals with the lowest order Bessel functions jo(x) and j±(x) are of 
great interest in applications related with the decays and photodetachment of atoms/ions 
in their ground states. The results of numerical computations of some integrals 
Ffc ; z;n(a, fi, 7), Bfcl n (a, (3, 7; V) and B^j. n (a, f3, 7; V) can be found in Tables I and II. For all 
integrals with Bessel functions shown in Table II we restricted to the accuracy ~ 1 ■ 10~ 15 . 
To determine the integrals B^J. n (a, /3, 7; V) and B^l n (a, (3, 7; V) to such an accuracy for 
V < 1 it was sufficient to use 30 terms in Eqs. (|T9|) and (120"]) . For V = 2 we used up to 75 
terms in these formulas. 
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The formulas for the three-body integrals Bf!^. n (a, /3, 7; V) with other spherical Bessel 
functions j^Vr^) can be obtained by using the same procedure. The result is 

ft 7; V) = V L g si^jiW., ft 7) (21) 

where k is integer number. In Eq. OTT) we have used the following formula for the 
Bessel function 

00 / -j\fc^2fe 

j L (z) = z L g 2kk[[2L - 2k - i)n (22) 

Dased on the use of the well-known 



13[ and analytical formulas for the 



Another approach for derivation of these formulas is 
recursion formulas for the spherical Bessel functions 
three-body integrals containing the jo(Vr 32 ) and ji(Vr 32 ) Bessel functions, Eqs.t TH?]) and 
(12 Op . The formulas for the integrals containing the spherical Bessel functions of different 
arguments, e.g., jx(Vr3i) and/or ji(Vr 2 i), are easily obtained from the expressions for the 
integrals with the iL(V A r 32 ) functions by applying a set of different a — > f3 — > 7 substitu- 
tions. Here we do not want to produce these formulas, since it is rather a technical problem 
which step-by-step repeates the procedure described above for deriving the formulas for 
the integrals with the joiVr^) and jiiVr^) Bessel functions. Instead we consider a more 
interesting and actual problem which is closely related with analytical and numerical compu- 
tation of other three-body integrals with modified Bessel functions Ko(r), Ki\(r) and K%2{r) 
functions. This problem is considered in the next Section. 



V. MATRIX ELEMENTS OF THE UEHLING POTENTIAL 



As is well known (see, e.g., 15], 16j) in the lowest order approximation the effect of 



vacuum polarisation between two interacting electric charges is described by the Uehling 



potential U(r) [17| . In [l8( we have derived the closed analytical formula for the Uehling 
potential. For atomic systems this formula is written in the following three-term form (in 
atomic units Ti = 1, m e 



3nr 



U(2br) 
(1 + 



l,e = l) 

2aQ 1 r f+™ 



2aQv^ b 2 r 2 



3tt 
K (2br 



exp(-2o; _1 ^r) (l + — A s 



br 
~6 



Ki x {2br) 



2i 
b 2 r 2 5 

— + 6 



1 



Ki 2 {2br) 



(23) 
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where the notation Q stands for the electric charge of the nucleus, b = a 1 and a = j- ~ jj= 
is the dimensionless fine-structure constant. Here and below % — is the reduced Planck 
constant (also called the Dirac constant), e is the electric charge of the positron and m e is 
the mass of the electron (= mass of the positron). In Eq.( l23l) Kq(o) is the modified Bessel 
function of zero order (see, e.g, [3]), i.e. 

■2k 



Z 

■ Z) 



POO ^ 

K (z) = exp(—zco8ht)dt = ^($(k + l)+]ii2-]ii~.., l , . , 
Jo k=0 2 2 "{k\y 

where i/j(k) is the Euler psi-function defined by Eq. (8.362) from 13] . The functions Kii(z) 
and Ki 2 (z) in Eq. f[2"3"j) are the recursive integrals of the K (z) = K (z) function, i.e. 

/*oo /*oo 

Kii(z) = / Ki (z)dz , and Ki n (z) = / Ki n ^i(z)dz , (24) 

J z J z 

where n > 1. 

The calculation of the matrix elements of the Uehling potential with the use of our three 
term formula, Eq. (l23l) . leads to the following three-body integrals: 

<(«,AT) = // jr^AK^M-^ - Pr 31 - 7^0****** (25) 

for p = 0,1,2 and two other similar integrals which contains Ki p (r 3 i) and Ki p (r 2 i). This 
problem can be solved approximately by using the known analytical formulas for the modified 
Bessel function Kq(z) and for the two lowest recursive integrals of this function. However, 
the overall accuracy of the final solution is not very high. There are some advanced methods 
which can be used to solve this problem, but at this moment we also trying to apply a few 
alternative approaches. One of these methods is based on the original integral representation 
for the Uehling potential, Eq.(]23~]). Below, we consider this procedure in details. 

First, consider the matrix elements of the regular Yukawa- type interparticle potential 
which is written in the form Vy(r) = Vq exp (~^ . The matrix element between the two 
exponential basis functions takes the following form 

V J J J r 31 r 2 i exp[-(a + fi)r 32 - f3r 31 - ir 2 i\dr m dr Z idr 2 i = V To-i-i(a + //, /3, 7) (26) 

in the case of the Vy{r 32 ) interaction. Note also that each of the integrals ro ; i ; i contains only 
four terms. This allows one to obtain the following formula in the case of the U 2 \ potential, 
which is the (21)-component of the total Uehling potential, Eq. f[23|) : 



C/ 2 i(2&£) =111 6XP<K ^ T21 K xp{-ar 32 - f3r 31 - ^r 21 )r 32 r 3l r 21 dr 32 dr 3l dr 21 
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(a + (3)(a + 7 + 26£)(/? + 7 + 26£) L (a + /?) 2 + (/? + 7 + 2&0(a + /?) + (2?) 

1 1 

+ 



(a + 7 + 2&0(a + /3) (P + 7 + 26^) (a + 7 + 2&£)- 
where a + /3>0,a + 7>0,/3 + 7>0 and £ > 0. The first factor 2 in the numerator is the 
Jacobian of the linear transformation from relative to perimetric coordinates. Analogous 
expressions can be obtained for the U{r 32 ) and t/(r 31 ) Yukawa-type potentials. In fact, such 
formulas can be derived from Eq.( l2"7|) simply by performing cyclic permutations of three 
parameters a, j3 and 7. 

The final formula for the matrix elements of the Uehling potential in an arbitrary Coulomb 
three-body system is written in the form 

2a r°°y „ „ „ i 1 \V?'-1 



37T Jl 



f™[qiq*U 21 (2b£) + qmU 31 (2bO + 9 2 g 3 ^32(2&0] (l + ja) 2 ^— # ' ( 28 ) 



where qi (i = 1, 2, 3) are the particle charges expressed in atomic units. The expressions 
for the Uij(2b£) terms are obtained from Eq. fl2T|) . By using the formulas, Eqs. fl27l) and 
(1281) . we have developed a number of effective numerical methods for accurate evaluation of 
the arising three-particle integrals. However, such methods are effective only for small and 
relatively small number(s) of basis functions used in calculations. Currently, our main goal 
is to derive the closed analytical expressions for all terms in Eq. (l28]) . since these formulas 
can be used at arbitrary number(s) of basis functions. 

VI. THREE-PARTICLE INTEGRALS WITH TWO BESSEL FUNCTIONS 

A number of actual problems in modern atomic physics lead to the three-particle integrals 
with two Bessel functions. For instance, the probability of formation of the negatively 
charged tritium ion (T~) during the (n, 3 He; 1 H, 3 H) nuclear reaction in the two-electron 
3 He-atom is reduced to the calculation of the following integral ('probability amplitude' for 



more details, see, e.g., [191 ]) 

Af = J J J $t- (r 32 , r 3 i, r 2 i)j (V t • r 32 )j (V t ■ r 3 i)^ He (^32, r 31 , r 2 i)r 3 2r3ir2idr 32 dr 31 dr 21 
1 

where V t is the speed of the tritium nucleus after the nuclear reaction in the 3 He atom. 



772 // / $t- (?"32, r 3 i, r 2 i)sin{V t ■ r 32 )sin(V t ■ r 31 )^ He (r 32 ,r 3 i, ^1)^1^32(^31(^21(29) 



To obtain the formula, Eq.(J29J), we have used the known fact from atomic physics that the 
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negatively charged hydrogen ion has only one bound 1 1 S'(L = 0)— state. Also, to derive 
Eq.f l29|) we restrict ourselves to the case when the incident 3 He atom was in its ground 
l l S{L = 0)-state. 

As mentioned in the second Section the wave functions of the ground l 1 5'(L = 0)— states 
in the two-electron H~ ion and He atom are usually approximated with the use of highly 
accurate variational expansion written in the relative/perimetric coordinates. The most 
advanced of such expansions is the exponential variational expansion in relative/perimetric 
coordinates which takes the following form (for the bound S(L = 0)— states in the three-body 
systems) : 

N 

Tp(r 32 , r 3 i, r 2i ) = ^ Qexp^a^ - A^3i - 7i r 2i) (30) 

i=l 

where Cj are the variational coefficients and N is the total number of terms in the trial 
function ip(r^2, T31, r?\) which is an accurate approximation of the actual wave function 
^( r 32i r 3i) r 2i)- The probability amplitude Aif is written as the double sum of the following 
three-particle integrals 

*g£<«. A 7) = / / / - fir 31 - 7 r 21 ) S mW • ^)sin(V t . r 31 )r M ^ m 

Theory of these integrals has not been developed in earlier studies. It was shown in 
that the integral, Eq.t lHTj) . is reduced to the following double sum (here we apply the Cauchy 
formula) 

B^olia, f3, 7; V) = fo ~ „x, Y C2 2 «i 1 rfe+2 M +i;«+2 K -2 At +i;n+i(a, P, l) (32) 
K =o \ AK + A )- fi=0 

where is the number of combinations from n by k (n > k) and i\z ;ra (a, b, c) is the basic 
three-particle integral defined above. 

The integral B^^a, /3, 7; V") belongs to the new class of three-particle integrals 
Bfcln(a, ft, 7; V). Such integrals contain the two Bessel functions jL{Vr 32 ) and ji{Vr 3 i). 
In this Section we consider one approach developed recently for computations of such inte- 
grals. First, by using the formula, Eq. (l2~2"|) . one finds the following formula for the product 
of the two Bessel functions j^ax) and jL 2 (fo/): 

00 (_-i\p P 2q 2qh2p~2q 2p-2q 

- g ^ E q (2Li + 2q + 1)!!(2L2 1 2p _ 2q + 1)H 03) 
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where p and q are both integer non- negative numbers and is the binomial coefficient. 
In our case we have a = b = V, x = r%2 and y = r^\. Therefore, the formula for the 
B^i-n" 2 \ a -> P\ 7i ^0 integrals take the form: 

B [LlL2) (o 8 -V- V) = V Ll+L2 V 1 ' V -5 x 

fc;Z;n 1 ' P ' 7 ' J ^ 2p P ! ^ (2L 1 + 2 ? + l)!!(2L 2 + 2p-2g+l)!! 

A 7) (34) 

This formula is appropriate in all applications where V < 10. For V < 1 the overall 
convergence rate of Eq. (|34l) for the B^^ 2 (a, 8, 7; V) integrals is fast, while for 1 < V < 2 
it is relatively fast and for 5 < V < 10 the convergence rate can be considered as moderate. 
In the cases when V > 15 — 20 one needs to develop some other methods of computation of 
the Bfcfin(a, (3, 7; V) integrals. Here we cannot discuss this interesting problem. 



VII. THREE-PARTICLE INTEGRALS FROM MORE COMPLICATED FUNC- 
TIONS 

By using the formulas derived above for three-particle integrals with one and two Bessel 
functions we can obtain some useful formulas for more complicated three-particle integrals. 
In reality, one finds a large number of three-particle integrals which can be approximated by 
the integrals with Bessel function(s). In this study we restrict ourselves to the consideration 
of the following integral: 

J(a,0,j;t) = J J f cos\]rl 2 - 2tr 32 • exp(-ar 32 - fir 31 - 7^21)^32^31^21 (35) 

Such integrals and their t— derivatives arise, e.g., in the problem of electron scattering on 
the electric dipole formed by the two heavy, positively charge particles. The formula for 
numerical evaluation of the J (a, 8, 7; t) integral is written in the form 

OO J.K 

J(a, p, 7; t) = £ -jBftl^a, /?, 7; 1) (36) 

All terms with k > 1 in this formula are determined directly with the use of formulas given 
above. Calculations of the term with k = contains a few additional troubles, since in 
this case k — 1 = — 1 and we need to define the integral S^.j. n (ot, /3, 7; 1). By using the 
formulas (10.1.11) and (10.1.12) from [14] one finds xj-i(x) = jo(x) — xji(x). Therefore, for 
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the •Sjtt-H. fl (a!, P, 7; 1) integral we have: 

si+iU «> Z 3 ' r, i) = A 7; i) - 32 M; >, & 7; 1) (37) 

Now, we can use the formula, Eq. fl36|) . to approximate the integral J (a, (3, 7; t) to arbitrarily 
high, in principle, numerical accuracy. There are many other useful for the integrals con- 
taining one and two Bessel functions and some functions of relative coordinates r 32 ,r 31 and 
r 21 . A number of such formulas will be considered in our next study. 

VIII. CONCLUSION 

We have considered the problems of analytical and numerical computation of the three- 
body (exponential) integrals which contain different Bessel functions. For a number of 
such integrals we have derived the closed analytical formulas and/or developed effective 
numerical methods. Our main interest is related to the integrals which contain the jo{krij) 
and ji(krij) Bessel functions, since such integrals play a central role in various problems on 
photodetachment of the ground S(L = 0)— states of different atomic and molecular systems. 
The formulas for the integrals with the spherical Bessel functions jiikrij) for L > 2 can 
be derived by applying the same procedure. The case of integrals with the modified Bessel 
functions, e.g., with the K n (brij) functions, is more difficult for analytical consideration. 
However, such integrals are needed in various problems, including derivation of the closed 
formulas for the matrix elements of the Uehling potential. 

Derivation of simple analytical formulas for the three-particle integrals with one and 
two Bessel functions has a great value in numerous applications to atomic and nuclear 
physics. On the other hand, it is very interesting mathematical problem, since three relative 
coordinates r 32 = r 2 , r 31 = r\ and r 12 always form a triangle. In fact, we are dealing with the 
new class of multiple integrals which have special form. Such integrals play a central role 
in various problems of three- and few-body physics. Note here that three-particle integrals 
with one and two Bessel functions are of paramount importance for atomic analysis of the 
species arising during the nuclear (n; £)— , (n; a) — and (n; p)— reactions in few electron atoms. 
Another interesting problem is calculation of the three-body integrals for the functions which 
are approximated by their series explicitly written in the spherical Bessel functions (and 
other Bessel functions). In earlier studies numerical evaluation of such integrlas was a very 
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difficult problem which had no effective and reliable solution. 

Appendix 

Addition theorem for the spherical Bessel functions 

Let us discuss here the statement known as the addition theorem for the spherical Bessel 
functions. This theorem plays a great role in the physics of three-body systems. On the 
other hand, it is an interest for the general theory of Bessel functions. First, consider the 
familiar Rayleigh expansion of a plane wave 

oo 

exp(ik ■ r 21 ) = £ i L (2L + l)j L (kr)P L (cosQ 21 ) (38) 

L=0 

where Jl{x) is the spherical Bessel function of scalar argument x, a - b designates the scalar 
product of the two vectors a and b, while Pl{v) is the Legendre polynomial. In Eq.(l38|) the 
notation G21 stands for the angle between the k and r 2 i vectors. 

Now, suppose we have the three-body system (123) and we need to re- write the Rayleigh 
expansion, Eq. f !38j) . in a different form which contains the plane waves of the two 'new' 
scalar products k ■ r 31 and k ■ r 32 . Such problems always arise in the three-body physics. 
It is clear that these two plane waves must uniformly be related to each other. The goal 
of this Appendix is to investigate such a relation in detail. Since the three particles in any 
three-body system always form a triangle, then we can write r 2 i = r 31 — r 32 . Now, one finds 
from Eq.flMD 

exp(ik ■ r 21 ) = exp(ik ■ r 3 i)exp(— tk ■ r 32 ) (39) 

By using the Rayleigh expansions twice for the right-hand side of the last equation we obtain 

00 00 00 

Y J i L {lL + l)j L {kr 21 )P L {cos® 21 ) = ^ Ll ^L 1 + l)j Ll {kr 31 )P Ll {cosQ 31 ) ■ £ i L * x 

L=0 Li=0 L 2 =0 

(2L 2 + l)j L2 (kr 32 )P L2 (-cosQ 32 ) (40) 

where the notations G 3 i and G 32 designate the angles between the k vector and r 3 i and r 32 
vectors, respectively. 

The right-hand side of the last equation can be transformed with the use of the Cauchy 
product defined by a discrete convolution as follows: 

00 L 

]T [£t L (-l) L - e (2£+ 1)(2L - 2£+l) 3e (kr 31 )j L . e (kr 32 )P e (cose 31 )P L ^(cose 32 )] (41) 
l=o e=o 
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It is crucially important here that each of the series in the Rayleigh expansion of a plane 

| I 

wave see, Eq. fl38|) and Eq.f l39|) ) converge absolutely [20]. From here and Eq. fHOj) we can 
write 

L 

(2L+l)j L (kr 21 )P L (cosQ 21 ) = (-l) L Y / (-lYW + V(2L-2e+l) 3£ (kr 31 )j L „ e (kr 32 ) x 

P e (cos0 31 )P L _ e (cosQ 32 ) (42) 

By multiplying the both sides of the last equation by Pl(cosQ 21 ) and performing the in- 
tegration over the 02i angle (or over all possible orientations of the k vector) we find the 
following expression 

L 

j L {kr 21 ) = (-l) L Y / (-lYW+l)(2L-2e+l) Je (kr 31 )j L ^(kr 32 ) x 

1=0 

j> P £ (cose 31 )P L ^(cosB 32 )P L (cose 21 )dk (43) 

This formula contains the spherical Bessel functions Jl(^^2i) written in terms of the spherical 
Bessel functions j n {kr 32 ) and j n {kr 3 i) of two other arguments, where the three vectors 
r 32 , r 3 i, r 2 i form a triangle, i.e. r 32 + r 2 i = r 3 i. As it follows from the definition of the cosQij 
angles we always have 02i + 631 + 632 = n. The formula, Eq.( 143l . is known as the addition 
theorem for the spherical Bessel functions. 
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TABLE I: Numerical values of some three-body integrals Tfc ; ; ;n (a, /?, 7). 



k I m a [i 7 r fc; ; ;n (a,/3,7) 7 r fe;Z;n (a, /3, 7) 

2 1 2.35 1.41 0.567 0.132484880489827E+00 -0.567 0.484535355001714E+01 

1 2 1 2.35 1.41 0.567 0.105479781157007E+00 -0.567 0.462617958966529E+01 

2 2 1 2.35 1.41 0.567 0.123759737118974E+00 -0.567 0.683620356276100E+01 

3 2 1 2.35 1.41 0.567 0.190628938378487E+00 -0.567 0.138242778966704E+02 

4 2 1 2.35 1.41 0.567 0.362095286177389E+00 -0.567 0.356816617385975E+02 

5 2 1 2.35 1.41 0.567 0.815657409095427E+00 -0.567 0.112342033402992E+03 

6 2 1 2.35 1.41 0.567 0.212162348108085E+01 -0.567 0.417926993577783E+03 

7 2 1 2.35 1.41 0.567 0.625059393550668E+01 -0.567 0.179435469496013E+04 

8 2 1 2.35 1.41 0.567 0.205551903374530E+02 -0.567 0.873301210942717E+04 

9 2 1 2.35 1.41 0.567 0.745934650018583E+02 -0.567 0.475056243580342E+05 



TABLE II: Numerical values of some three-body integrals B^.j.^a, f3, 7; V) and B k . ln (a, (3, 7; V). 

klm a f) 7 V B^a, f3, T ,V) B^ja, /3, r ,V) 

3 2 1 2.35 1.41 0.567 0.25 0.18233241643012516E+00 0.290930992106451E-01 

5 2 1 2.35 1.41 0.567 0.25 0.75291471135429875E+00 0.166432412830887E+00 

3 2 1 2.35 1.41 0.567 0.50 0.15968050735256670E+00 0.522255954684081E-01 

5 2 1 2.35 1.41 0.567 0.50 0.59041249572520414E+00 0.278021233893212E+00 

3 2 1 2.35 1.41 0.567 1.00 0.94868174980045456E-01 0.691516883096556E-01 

5 2 1 2.35 1.41 0.567 1.00 0.20605506256710767E+00 0.274928833359198E+00 

3 2 1 2.35 1.41 0.567 1.50 0.40374337963233781E-01 0.554457473644749E-01 

5 2 1 2.35 1.41 0.567 1.50 0.20605506256710767E+00 0.274928833359198E+00 

3 2 1 2.35 1.41 0.567 2.00 0.11173049407361310E-01 0.340384106321226E-01 

5 2 1 2.35 1.41 0.567 2.00 -0.35522376544132919E-01 0.316198754574614E-01 
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